Pre-visualization of the normalized count data before any differential analysis

Load¶

Libraries and functions¶

Warning message in is.na(x[[i]]):
“is.na() applied to non-(list or vector) of type 'environment'”Warning message in rsqlite_fetch(res@ptr, n = n):
“Don't need to call dbFetch() for statements, only for queries”

==========================================================================
*
*  Package WGCNA 1.63 loaded.
*
*    Important note: It appears that your system supports multi-threading,
*    but it is not enabled within WGCNA in R. 
*    To allow multi-threading within WGCNA with all available cores, use 
*
*          allowWGCNAThreads()
*
*    within R. Use disableWGCNAThreads() to disable threading if necessary.
*    Alternatively, set the following environment variable on your system:
*
*          ALLOW_WGCNA_THREADS=<number_of_processors>
*
*    for example 
*
*          ALLOW_WGCNA_THREADS=4
*
*    To set the environment variable in linux bash shell, type 
*
*           export ALLOW_WGCNA_THREADS=4
*
*     before running R. Other operating systems or shells will
*     have a similar command to achieve the same aim.
*
==========================================================================


Allowing multi-threading with up to 4 threads.
[1] "preparing gene to GO mapping data..."
[1] "preparing IC data..."
[1] "preparing gene to GO mapping data..."
[1] "preparing IC data..."
[1] "preparing gene to GO mapping data..."
[1] "preparing IC data..."

Data¶

Sample clustering¶

With all genes¶

Without X/Y genes¶

PCA¶

Updating HTML index of packages in '.Library'
Making 'packages.html' ... done

'data.frame':	150 obs. of  5 variables:
 $ Sepal.Length: num  5.1 4.9 4.7 4.6 5 5.4 4.6 5 4.4 4.9 ...
 $ Sepal.Width : num  3.5 3 3.2 3.1 3.6 3.9 3.4 3.4 2.9 3.1 ...
 $ Petal.Length: num  1.4 1.4 1.3 1.5 1.4 1.7 1.4 1.5 1.4 1.5 ...
 $ Petal.Width : num  0.2 0.2 0.2 0.2 0.2 0.4 0.3 0.2 0.2 0.1 ...
 $ Species     : Factor w/ 3 levels "setosa","versicolor",..: 1 1 1 1 1 1 1 1 1 1 ...

Gene co-expression network analysis¶

Build network¶

pickSoftThreshold: will use block size 4378.
 pickSoftThreshold: calculating connectivity for given powers...
   ..working on genes 1 through 4378 of 10219
   ..working on genes 4379 through 8756 of 10219
   ..working on genes 8757 through 10219 of 10219
   Power SFT.R.sq  slope truncated.R.sq  mean.k. median.k.  max.k.
1      2    0.144 -0.617          0.901 644.0000  6.20e+02 1390.00
2      4    0.691 -1.460          0.964 114.0000  9.21e+01  438.00
3      6    0.783 -1.840          0.970  31.1000  1.94e+01  198.00
4      8    0.790 -2.020          0.962  11.0000  5.03e+00  106.00
5     10    0.832 -1.990          0.985   4.6100  1.49e+00   63.10
6     12    0.851 -1.960          0.988   2.1900  4.90e-01   40.10
7     14    0.876 -1.890          0.988   1.1500  1.72e-01   26.70
8     16    0.886 -1.820          0.971   0.6430  6.43e-02   18.40
9     18    0.899 -1.780          0.970   0.3830  2.50e-02   13.60
10    20    0.882 -1.770          0.946   0.2400  1.02e-02   10.40
11    22    0.917 -1.720          0.984   0.1560  4.29e-03    8.11
12    24    0.902 -1.700          0.962   0.1050  1.86e-03    6.46
13    26    0.408 -2.230          0.416   0.0724  8.32e-04    5.23
14    28    0.416 -2.190          0.427   0.0514  3.76e-04    4.28
15    30    0.931 -1.630          0.984   0.0372  1.74e-04    3.55

Parameters for network construction

Soft-thresholding power parameter: 16
Minimum module size: 30
Height for module merging: 0.2 (below)
Biweight Midcorrelation

ERROR while rich displaying an object: Error in sprintf(wrap, header, body): 'fmt' length exceeds maximal format length 8192

Traceback:
1. FUN(X[[i]], ...)
2. tryCatch(withCallingHandlers({
 .     rpr <- mime2repr[[mime]](obj)
 .     if (is.null(rpr)) 
 .         return(NULL)
 .     prepare_content(is.raw(rpr), rpr)
 . }, error = error_handler), error = outer_handler)
3. tryCatchList(expr, classes, parentenv, handlers)
4. tryCatchOne(expr, names, parentenv, handlers[[1L]])
5. doTryCatch(return(expr), name, parentenv, handler)
6. withCallingHandlers({
 .     rpr <- mime2repr[[mime]](obj)
 .     if (is.null(rpr)) 
 .         return(NULL)
 .     prepare_content(is.raw(rpr), rpr)
 . }, error = error_handler)
7. mime2repr[[mime]](obj)
8. repr_latex.matrix(obj)
9. repr_matrix_generic(obj, sprintf("\\begin{tabular}{%s}\n%%s%%s\\end{tabular}\n", 
 .     cols), "%s\\\\\n\\hline\n", "  &", " %s &", "%s", "\t%s\\\\\n", 
 .     "%s &", " %s &", escape_fun = latex_escape_vec, ...)
10. sprintf(wrap, header, body)

 Calculating module eigengenes block-wise from all genes
   Flagging genes and samples with too many missing values...
    ..step 1
 ....pre-clustering genes to determine blocks..
   Projective K-means:
   ..k-means clustering..
   ..merging smaller clusters...
Block sizes:
gBlocks
   1    2    3 
4997 4015 1207 
 ..Working on block 1 .
    TOM calculation: adjacency..
    ..will use 4 parallel threads.
     Fraction of slow calculations: 0.000000
    ..connectivity..
    ..matrix multiplication (system BLAS)..
    ..normalization..
    ..done.
   ..saving TOM for block 1 into file norm_genes_TOM-block.1.RData
 ....clustering..
 ....detecting modules..
 ....calculating module eigengenes..
 ....checking kME in modules..
 ..Working on block 2 .
    TOM calculation: adjacency..
    ..will use 4 parallel threads.
     Fraction of slow calculations: 0.000000
    ..connectivity..
    ..matrix multiplication (system BLAS)..
    ..normalization..
    ..done.
   ..saving TOM for block 2 into file norm_genes_TOM-block.2.RData
 ....clustering..
 ....detecting modules..
 ....calculating module eigengenes..
 ....checking kME in modules..
 ..Working on block 3 .
    TOM calculation: adjacency..
    ..will use 4 parallel threads.
     Fraction of slow calculations: 0.000000
    ..connectivity..
    ..matrix multiplication (system BLAS)..
    ..normalization..
    ..done.
   ..saving TOM for block 3 into file norm_genes_TOM-block.3.RData
 ....clustering..
 ....detecting modules..
 No modules detected in block 3
  ..reassigning 1 genes from module 1 to modules with higher KME.
  ..reassigning 1 genes from module 9 to modules with higher KME.
  ..reassigning 3 genes from module 10 to modules with higher KME.
  ..reassigning 3 genes from module 11 to modules with higher KME.
  ..reassigning 1 genes from module 14 to modules with higher KME.
 ..merging modules that are too close..
     mergeCloseModules: Merging modules whose distance is less than 0.2
       Calculating new MEs...

   0    1    2    3    4    5    6    7    8    9   10   11 
8786  446  260  191   95   93   86   62   56   53   49   42

11 modules with between 42 and 446 genes
8786 genes not assigned to a module

Network visualization¶

TOM calculation: adjacency..
..will use 4 parallel threads.
 Fraction of slow calculations: 0.000000
..connectivity..
..matrix multiplication (system BLAS)..
..normalization..
..done.

Module analysis¶

[1] 1

Warning message in grep(txdbPattern, installedPackages):
“argument 'pattern' has length > 1 and only the first element will be used”Warning message in pcls(G):
“initial point very close to some inequality constraints”Fetching GO annotations...
For 2918 genes, we could not find any categories. These genes will be excluded.
To force their use, please run with use_genes_without_cat=TRUE (see documentation).
This was the default behavior for version 1.15.1 and earlier.
Calculating the p-values...
'select()' returned 1:1 mapping between keys and columns

        category over_represented_pvalue under_represented_pvalue numDEInCat
3121  GO:0006325            1.535024e-17                        1         59
5998  GO:0016569            1.552323e-14                        1         48
4834  GO:0010467            2.180602e-13                        1        177
5999  GO:0016570            8.098099e-13                        1         40
9705  GO:0034641            9.689717e-13                        1        194
14200 GO:0051171            1.033361e-12                        1        156
      numInCat                                              term ontology
3121       683                            chromatin organization       BP
5998       489                            chromatin modification       BP
4834      4726                                   gene expression       BP
5999       384                              histone modification       BP
9705      5508      cellular nitrogen compound metabolic process       BP
14200     3783 regulation of nitrogen compound metabolic process       BP
[1] 2

Warning message in grep(txdbPattern, installedPackages):
“argument 'pattern' has length > 1 and only the first element will be used”Warning message in pcls(G):
“initial point very close to some inequality constraints”Fetching GO annotations...
For 2918 genes, we could not find any categories. These genes will be excluded.
To force their use, please run with use_genes_without_cat=TRUE (see documentation).
This was the default behavior for version 1.15.1 and earlier.
Calculating the p-values...
'select()' returned 1:1 mapping between keys and columns

        category over_represented_pvalue under_represented_pvalue numDEInCat
9908  GO:0034976            3.637326e-08                1.0000000         16
3468  GO:0006810            1.089150e-07                1.0000000         94
5357  GO:0015031            1.267559e-07                1.0000000         50
14240 GO:0051234            3.386036e-07                1.0000000         95
12359 GO:0045184            6.817248e-07                0.9999997         52
13425 GO:0048193            9.479241e-07                0.9999998         14
      numInCat                                     term ontology
9908       204 response to endoplasmic reticulum stress       BP
3468      4082                                transport       BP
5357      1649                        protein transport       BP
14240     4231            establishment of localization       BP
12359     1809    establishment of protein localization       BP
13425      186                  Golgi vesicle transport       BP
[1] 3

Warning message in grep(txdbPattern, installedPackages):
“argument 'pattern' has length > 1 and only the first element will be used”Fetching GO annotations...
For 2918 genes, we could not find any categories. These genes will be excluded.
To force their use, please run with use_genes_without_cat=TRUE (see documentation).
This was the default behavior for version 1.15.1 and earlier.
Calculating the p-values...
'select()' returned 1:1 mapping between keys and columns

        category over_represented_pvalue under_represented_pvalue numDEInCat
3468  GO:0006810            2.086915e-24                        1         82
14240 GO:0051234            2.515869e-24                        1         83
14206 GO:0051179            6.256614e-24                        1         92
14955 GO:0055114            9.097327e-22                        1         43
12773 GO:0046034            6.631370e-20                        1         24
4467  GO:0009167            3.481026e-19                        1         24
      numInCat                                                  term ontology
3468      4082                                             transport       BP
14240     4231                         establishment of localization       BP
14206     5342                                          localization       BP
14955      883                           oxidation-reduction process       BP
12773      181                                 ATP metabolic process       BP
4467       206 purine ribonucleoside monophosphate metabolic process       BP
[1] 4

Warning message in grep(txdbPattern, installedPackages):
“argument 'pattern' has length > 1 and only the first element will be used”Warning message in pcls(G):
“initial point very close to some inequality constraints”Fetching GO annotations...
For 2918 genes, we could not find any categories. These genes will be excluded.
To force their use, please run with use_genes_without_cat=TRUE (see documentation).
This was the default behavior for version 1.15.1 and earlier.
Calculating the p-values...
'select()' returned 1:1 mapping between keys and columns

[1] category                 over_represented_pvalue  under_represented_pvalue
[4] numDEInCat               numInCat                 term                    
[7] ontology                
<0 rows> (or 0-length row.names)
[1] 5

Warning message in grep(txdbPattern, installedPackages):
“argument 'pattern' has length > 1 and only the first element will be used”Warning message in pcls(G):
“initial point very close to some inequality constraints”Fetching GO annotations...
For 2918 genes, we could not find any categories. These genes will be excluded.
To force their use, please run with use_genes_without_cat=TRUE (see documentation).
This was the default behavior for version 1.15.1 and earlier.
Calculating the p-values...
'select()' returned 1:1 mapping between keys and columns

[1] "1: 12 - 2.69058295964126%"
[1] "2: 14 - 5.38461538461539%"
[1] "3: 6 - 3.1413612565445%"
[1] "4: 1 - 1.05263157894737%"
[1] "5: 1 - 1.0752688172043%"
[1] "6: 4 - 4.65116279069767%"
[1] "7: 2 - 3.2258064516129%"
[1] "8: 2 - 3.57142857142857%"
[1] "9: 2 - 3.77358490566038%"
[1] "10: 0 - 0%"
[1] "11: 0 - 0%"

Functional gene network¶

Warning message in scan(file = file, what = what, sep = sep, quote = quote, dec = dec, :
“EOF within quoted string”Warning message in scan(file = file, what = what, sep = sep, quote = quote, dec = dec, :
“number of items read is not a multiple of the number of columns”Warning message in scan(file = file, what = what, sep = sep, quote = quote, dec = dec, :
“EOF within quoted string”Warning message in scan(file = file, what = what, sep = sep, quote = quote, dec = dec, :
“number of items read is not a multiple of the number of columns”Warning message in scan(file = file, what = what, sep = sep, quote = quote, dec = dec, :
“EOF within quoted string”Warning message in scan(file = file, what = what, sep = sep, quote = quote, dec = dec, :
“number of items read is not a multiple of the number of columns”

Tree with all GOs at "../results/dge/global_GO_network.png"

Median of gene counts in the modules over ages¶

Module 1: chromatine/chromosome organization/RNA metabolic process

Module 2: response to endoplasmic reticulum stress/transport

Module 3: metabolic process (ATP, ribonucleoside)

Module 4

Module 5

Module 6: metabolic process (primary, cellular)

Module 7: immune system

Module 8: translation/rRNA

Module 9

Module 10: organelle

Module 11: localization

	GF_104w_F_1_2	GF_104w_F_2_2	GF_104w_F_3_2	GF_104w_M_1_2	GF_104w_M_2_2	GF_52w_M_1_2	GF_52w_M_2_2	GF_52w_M_3_2	GF_52w_M_4_2	GF_8w_M_1_2	⋯	SPF_104w_M_14_2	SPF_8w_F_1_2	SPF_8w_F_3_2	SPF_8w_F_4_2	SPF_8w_F_5_2	GF_8w_F_1_2	GF_8w_F_2_2	GF_8w_F_3_2	GF_8w_F_4_2	GF_8w_F_5_2
0610005C13Rik	0.0000	0.00000	0.0000	1.768593	1.259244	0.8418979	4.942785	5.603216	7.231584	8.47277	⋯	5.28524	3.500355	2.112002	2.602398	0.00000	18.24998	0.00000	2.025447	5.268672	0.00000
0610007P14Rik	868.4015	749.64927	910.3535	886.064952	992.913601	508.5063502	530.855142	514.375197	738.344705	603.98749	⋯	692.36648	763.077332	744.480801	908.237050	957.79030	898.54336	682.26034	612.697796	593.252445	809.27852
0610009B22Rik	1021.3276	1169.19730	1180.9991	1224.750457	1558.313990	1038.9020466	876.850113	1016.423320	877.914273	948.95029	⋯	1206.79653	1268.295198	1164.769254	1295.994415	1446.89758	1040.24912	1020.67017	1142.352254	959.952003	1094.00091
0610009L18Rik	28.2185	36.20465	18.2253	20.338816	38.406931	31.9921214	26.691041	22.412863	4.338950	61.73018	⋯	45.80542	42.004257	45.408049	102.794738	72.02239	48.30878	19.58642	11.139960	45.310578	24.36813
0610009O20Rik	217.5555	278.98879	266.0893	226.379868	271.996624	299.7156634	349.949200	395.587025	309.511786	360.69794	⋯	178.81730	274.194454	252.384272	227.709862	184.89330	231.88216	158.86764	271.409933	304.529230	269.33199
0610010B08Rik	0.0000	0.00000	0.0000	0.000000	0.000000	0.0000000	0.000000	0.000000	0.000000	0.00000	⋯	0.00000	0.000000	0.000000	0.000000	0.00000	0.00000	0.00000	0.000000	0.000000	0.00000

	0610007P14Rik	0610009B22Rik	0610009L18Rik	0610009O20Rik	0610010F05Rik	0610010K14Rik	Mettl26	0610012G03Rik	0610030E20Rik	0610037L13Rik	⋯	Zw10	Zwint	Zxda	Zxdb	Zxdc	Zyg11b	Zyx	Zzef1	Zzz3	Hikeshi
GF_104w_F_1_2	868.4015	1021.328	28.21850	217.5555	411.4439	1351.757	283.0953	693.6289	899.3508	449.6754	⋯	618.0761	1471.913	297.6596	130.1692	923.9282	1333.552	2864.633	396.8795	752.7967	668.1412
GF_104w_F_2_2	749.6493	1169.197	36.20465	278.9888	461.0769	1516.336	207.6443	686.8236	759.2329	526.0323	⋯	797.5672	1539.763	243.8490	198.0607	981.7850	1340.637	3358.514	364.1762	597.3768	703.8610
GF_104w_F_3_2	910.3535	1180.999	18.22530	266.0893	288.8709	1852.601	222.3486	868.4353	540.3800	444.6972	⋯	731.7456	1375.099	153.0925	171.3178	922.1999	1181.910	3058.205	209.5909	644.2642	770.9300
GF_104w_M_1_2	886.0650	1224.750	20.33882	226.3799	278.5534	1551.056	190.1237	823.2799	525.2720	607.5116	⋯	466.9085	1752.675	172.4378	344.8756	683.5611	1295.494	2789.955	261.7517	944.4285	687.9826
GF_104w_M_2_2	992.9136	1558.314	38.40693	271.9966	304.1073	1911.532	243.0340	778.8422	707.6949	586.8075	⋯	878.3224	1957.494	147.9611	249.3302	1003.6172	1205.726	2819.446	200.8494	785.1384	583.6594
GF_52w_M_1_2	508.5064	1038.902	31.99212	299.7157	392.3244	1173.606	182.6919	589.3286	1047.3210	621.3207	⋯	768.6528	1397.551	287.9291	275.3006	924.4039	1674.535	2206.614	749.2892	1163.5029	493.3522

	0610030E20Rik	1600020E01Rik	Zfp974	1700094D03Rik	1700112E06Rik	1810026B05Rik	⋯	Zmym6	Znhit2	Znrf1	Zscan26	Zswim8
0610030E20Rik	1	1	0	0	0	1	⋯	0	0	0	1	0
1600020E01Rik	1	1	0	1	0	1	⋯	1	0	1	1	1
Zfp974	0	0	1	0	0	0	⋯	0	1	0	0	0
1700094D03Rik	0	1	0	1	0	1	⋯	1	0	0	1	0
1700112E06Rik	0	0	0	0	1	0	⋯	0	0	0	0	0
1810026B05Rik	1	1	0	1	0	1	⋯	1	0	1	1	1

	GO:0050896	GO:0048731	GO:0048856	GO:0044707	GO:0032502	GO:0007275	⋯
0610009O20RIK	1	0	0	0	0	0	⋯
FMC1	0	0	0	0	0	0	⋯
1190002N15RIK	1	1	1	1	1	1	⋯
1600002K03RIK	0	0	0	0	0	0	⋯
1600014C10RIK	0	0	0	0	0	0	⋯
ZFP974	0	0	0	0	0	0	⋯

Sepal.Length	Sepal.Width	Petal.Length	Petal.Width	Species
5.1	3.5	1.4	0.2	setosa
4.9	3.0	1.4	0.2	setosa
4.7	3.2	1.3	0.2	setosa
4.6	3.1	1.5	0.2	setosa
5.0	3.6	1.4	0.2	setosa
5.4	3.9	1.7	0.4	setosa